Background

With this script, I want to be able to input which models to run and what type of output I want. The rest should be automated and pull from the sorucefiles that I have built my models from. Ideally, I want this to be the file that Tom can use in future, zipped up with all the dependencies so it is very easy to extract data.

Set working directories & Load sourcefiles

If this is running from the shared Rproject, this should all be automated.

Comparing model results

Specify what data you are interested in

  • models: Here you can input which models you want to compare with each other, make sure to put them between "" and don’t include , or . in the model numbers. E.g. use 11, 12 and 332. If you are just interested in 1 model, this is fine as well. Models x.1 are the non-hoarding models. x.2 are the left-over hoarders. models x.3.1 are direct hoarders that have hoarding as the ‘top’ behaviour (above highest threshold). Models x.3.2 are direct hoarders that have resting as the ‘top’ behaviour (above the highest threshold). As a default, this is set to all models, if you want to change this you need to run the code from scratch, which will take more time.

  • int_var: Here you specify your variable of interest. These are only relevant for the graphs that are split per enviornment, where you will visualise the chosen variable across the 12 environments

    • eat = proportion of alive birds that are eating at any given timestep
    • eat_hoard = proportion of alive birds that are doing the ‘eat and hoard leftovers’ behaviour at any given timestep. Note that this is only an option for X.2 models.
    • forage = proportion of alive birds that goes out to forage. Foraging happens before the ‘eating’ behaviour, before the ‘direct hoarding’ behaviour and before ‘eat_haording’ behaviour.
    • dir_hoard= the proportion of alive birds that does the ‘direct hoarding’ behaviour at any given timestep. Note that only X.3.Y models can do this.
    • alive = The percentage of birds that is alive at any given timestep
    • caches = The number of caches birds have at any give timestep. Note that non-hoarding species do get a starting level of cahches, they just cannot access them.
    • find_food = The number of food items that are found in each step
    • fat_res = The fat reserve of alive birds in each timestep
    • stom_con= The stomach content of alive birds in each timestep
    • fat_loss_r = The fat loss rate of alive birds in each timestep
    • mass= The mass of alive birds in each timestep
    • predation = Average of successful predation attempts in each timestep
    • rest = The proportion of alive birds that is resting in each timestep
    • retrieve = The proportion of alive birds that is retrieving in each timestep. Note this is only possible for hoarding models (X.2 and X.3.Y)
    • sleep= The proportion of alive birds that is sleeping in any given timestep
    • cache_decay = The average of decayed caches in every given timestep
    • temp = Temperature in degrees Celcius
  • num_birds: specify the number of individual birds you want to run each model for (default = 1000). Only change this if you want to run the models from scratch. You also need to change which code chunks are running in that case.

  • daylight_h: specify how many daylight hours you want (default = 8). Only change this if you want to run the models from scratch. You also need to change which code chunks are running in that case.

  • days: specify how many days you want in your simulation (default = 30). Only change this if you want to run the models from scratch. You also need to change which code chunks are running in that case.

Repeat current settings

The current run includes the following settings:

  • Models to be evaluated = 11, 21, 31, 12, 22, 32, 131, 231, 331, 132, 232
  • Variable of interest (if we are looking at split environments) = stom_con
  • Number of birds = 1000
  • Daylight hours = 8
  • Number of days in simulation = 30

Then follows the code that loops through each of the models specified and runs them for 30 days and 1000 individuals. In future, I’ll make this flexible at the top as well, but for now it is hard coded. The code is not printed here, but you can look at it in the Rmarkdown file. This is followed by code for the visualization. What is visualized depends on the above input variables.

setwd("C:/Local_R/BiPHD-ABM/")
load( "results.Rda" )

Survival across the selected models

The following graph shows how the different models survive across the environments. Survival is hte proportion of birds alive at any given timestep.

Mean variable values across all environments

This graph shows the average of each ofthe variables on each timestep. These averages are taken across all environments.

Variable of interest, split per enviornment

These graphs show the chosen variable of interest (see above) for each of the environments. The averages are taken across all alive birds in the simulation

Daily patterns

With this code I want to do the following:

  • Aggregate over daily patterns
  • Omit the first 3 days –> This does cause some models to not show up, or show up half way. The number of datapoints in each aggregation will also differ due to survival. We need to consider how we want to deal with that.
  • Compare the pattern of fat-loss, SC and FLR
plot<-ggplot(all_vars_day_pat_tot, aes(x=timestep_within_day, y=mean_new))+
        geom_line(aes( col=mod_id, linetype=line_id), size=1)+
        facet_wrap(.~id, scales='free_y',nrow=6)+
        coord_cartesian(xlim=c(0, 25))+
        ggtitle(label="All variables - Daily patterns - Across all environments")+
        labs(y=paste(int_var), x="Timestep within the Day", col="Model ID")+
        theme(plot.title=element_text(size=25, face='bold'), 
              axis.title.x = element_text(size=20, face='bold' ), 
              axis.title.y=element_text(size=20, face='bold'), 
              strip.text=element_text(size=20), 
              legend.text=element_text(size=20), 
              legend.title = element_text(size=20, face='bold'))+
  scale_color_manual(values=colours)


p<-ggplotly(plot)
p